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SUMMARY 


The .tetos of three tesks related to reaction rate. In 
turbulent, reactin, tloos are revieoed. The work dealing with ^ 
assumed - pdf modeling of reaction rates is 

ting a three-variable pdf en®>loying a "most-likely P ’ 

studiea related to appropriate chemical kinetic nechaniane treat- 
ing hydrogen-air co*ustion have been initiated, their curren 
atatu. la reviewed. Perfeotly-atirred reactor modeling of 
fl«»e-.tabilizing recirculation regiona waa uaed to Inveatigate 
the .table flame region, for allane, hydrogen, metnane an 
propane, and for certain mUturea theraof. It ia concluded t^t 
in general, ailane c«, be counted upon to atabilite fl«aea only 
When the overall fuel-air .guivalence ratio i. “ 

greater than unity. For lean flamea, ail«ie may 
the flame. Other factor, favoring .table flame, are high initial 

temoerature and system pressure. 
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1. Introduction 


1.1 Overview 

The prediction of turtclent reeding flows U . mejor concern of 
current combustion reseerch. Such modeling is required in order to enhence 
the understanding of the phenomena involved and to design and evaluate the per- 
formance of combustion devices. A principle element to be derived from the 
modeling of these flows is an expression for the reaction rates of the various 
species involved in any particular combustion process under consideration. 

Currently, several approaches tor the determination of the properties 
of turbulent, reacting flows exist. Ill Of the present approaches, the method of 
utilizing an assumed probability density function (pdf) for temperature and species 
concentrations is selected for use in this study. The motivation for selecting 

is its relative computational simplicity and basis in the probabilistic 


this approach is 
nature of turbulence. 


1.2 Objective 


This work examines the effects of temperature and species concentration 
fiuctua.ions on reaction rates in turbulent reacting flows by means of the assumed 


pdf approach. The following items are 


the subject of the present study. 
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1 . The most-likely joint pdf of reference [2] to describe effects of 
temperature and species concentration fluctuations on the reaction rate. 

2. The extension of the most-likely joint pdf of reference [2] to three 
variables. This pdf is also used to describe temperature and species concentration 
fluctuations on the reaction rate. 

Section 2 presents the general theory of pdf's and their use in calculating 
the mean turbulent reaction rate constant. Section 4 presents the results of 
these calculations. Section 5 presents the conclusions and direction for further 

study. 


OF POOR ^ 
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2. General Theory 


2.1 Continuous Random Variables and Probability Density Functions 

A variable x may be thought of as a continuous random variable if it 
can assume all values in some interval, where the endpoints of the interval may 
be plus and minus infinity. In a rigorous mathematical sense, x is said to be 
a continuous random variable if there exists a function p(x), which satisfies 
the following conditions: [31 


p(x) ^ o 


( 1 ) 


/* p(x)dx « 1 
- 00 ^ 

It should be noted that p(x) by itself does not represent a probability, but rather 
the area under the curve from a to b represents the probability that the continuous 

r 

random variable x will have a value somewhere on the interval. Thus the inter- 
pretation of equation (lb) is that the probability that x has a value somewhere over 
the entire range of values is unity. If the pdf is defined on a finite closed interval, 
[a,b], the probability of x having a value outside this interval must of course be 
zero. 


The expected, or mean value, of a continuous random variable x is 
defined as: [3] 

= E(x) soo^* X C2) 
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The variance of this same variable is defined as: [31 

\oJ ( * “ ^ 

m 

The square root of the variance is the standard deviation. 

The concept of the mean value of a continuous random variable can be 
extended to functions of a continuous random variable . For example , let h(x) 
be such a function. Then, the mean value of h(x) is: [3] 

h(x) = . oj/ h(x) p(x) dx (4) 

Probability density functions can be defined for more than one continuous 
random variable. These are sometimes termed joint, multivariable or multi- 
dimensional pdf's. For a joint pdf of two continuous random variables, the 
probability is represented by the volume under the surface described by the 

t 

pdf. For such a pdf, p(x,y), the conditions corresponding to equation (1) are: 
p(x,y) ^ o 

» (5) 

/ /“ p(x,y)dxdy = 1 

-00 *<» 

The coTicepts of ‘’moments" of a one-dimensional continuous random 
variable may be extended to multi-dimensional continuous random variables. For 
example, the joint moment about the origin of a two-dimensional continuous random 
variable, (x,y), is expressed as: 13] 

^ xy xyp (x,y) dxdy (6) 

By comparison, equation (6) is seen to be an extension of equation (2). 




The joint moment 
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about the mean of a two-dimensional continuous random 


variable, (x,y)t is expressed as: [3] 

V ?•/*!•/* (X - ECx) ) (y - E (y) ) p(x,y) dxdy (7) 

u is also called the covariance of x and y. 

xy 

By comparison, equation (7) is seen to be an extension of equation (3). 

If a joint pdf is specified, it is possible to examine the distribution of any 
one of the continuous random variables by consideration of its marginal pdf. 

For a two-dimensional pdf. p(x.y). the marginal pdf of x. for example is given 

by: 13) 

hCx) = ^„/*pCx,y)dy 

This marginal distribution of x may be thought of as the distribution of x, with 
the simultaneous behavior of the other variables suppressed. In other words, . 

only the behavior of x is being examined. 

Using the concept of a marginal pdi, the moment about the origin ot any of 
the one-dimensional components ot a multi-dimensional continuous random variable 
may be expressed. For a two-dimensional continuous random variable, the moment 

about the origin of x, for example, is given by: 

y(x) =.«/“«/* xp(x,y)dydx ^9) 

Examination of equations (8) and (9) reveals that the moment about the origin of x 

is expressedias the integral from - • to +» , with respect to x. of the product 
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of X and its marginal pdf. 

In the case of joint pdf's and functions of more than one continuous random 
variable, a similar expression to equation (4) can be written. For a function 
dependent of two continuous random variables, say g(x,y), the mean value is 
written as: [3] 

g(x,y) «.^/?«/*g(x»y)p(x,y)dxdy (9) 

The correlation coefficient, P xy» is a parameter defined for the two- 
dimensional continuous random variable (x,y) as: 

Pxy = Ef [x-E(x)l [y-E(y)3 } (10) 

V (x) V (y) 

The numerator of the correlation coefficient is the covariance, shown in equation 
(7). 

It is demonstrated in Reference I3] that the value of the correlation 
coefficient is always between -1 and +1 . 

If s ^ 1 , X and y are linear function of each other. A negative 

value of P xy indicates one variable increases while the other decreases. A 
positive value indicates that both variables increase or decrease simultaneously. 

A value of P close to zero only indicates the absence of a linear relationship between 
X and y. It does not preclude the possibility of some nonlinear relationship. 

2.2 Arrhenius Reaction Rate Constant and Reaction Rate Expressions 

For the purpose of this study, the following one-step, irreversible 


8 
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reaction is considered: 
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F + V O **■? 

F,0,P « fuel, oxidizer, and product species, respectively 
V m stoichiometric coefficient (molar) 


(12) 


The disappearance of F and 0 are related by: 


"<lCp “1 (13) 

dt V dt 

where: Cp , Cq * concentrations of fuel and oxidizer, respectively 

The negative signs in equation (13) indicate that the concentrations of F and O 
decrease as the reaction proceeds. The reaction rate of fuel is expressed as: 


*F 


~^^F = k(T)Cp''CQ*" 

dt 


(14) 


' where k(T) * Arrhenius reaction rale constant 
T » temperature , 

n,M ® constants dependent upon the particular reaction 

The Arrhenius reaction rate constant is expressed as: 

k(T) » AT®exp(-T^/T) ^^5) 

where A,B = constants dependent upon the particular reaction 

T = activation temperature of the particular reaction 
A 

Equation (15) can be modified by the introduction of the dimensionless temperature: 


t 


T 


max 



( 16 ) 
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Th. «n.p.r..ur.s .nd T„,„ .r. detlr.«i «. «h.. the dimensionles. ,.n.p.r..ur. 

, c.n onl, .ssum. ..lu« vdthln Ih. lnt.rv.1 (O.ll. Thus. T.,„ ».y b. <h. low«t 
.emeratur. of .h. unr«c«d »na,l.u.n.s and T„„ ma, b. th. «,«llibri«m cmbua.loa 

temperature of th. reaction. Sub.tliutlon of (16) into ( 15) yields: 

k(0 - AOtjt * kj)® exp l-V“'l' * V ' - 

Where kj » 

^2 * *^min 

Equation (U) is also modified by the introduction of the following dimensionless 


concentrations: 



( 18 ) 


In these expressions, "max" and "min" denote maximum and minimum values, 
respectively. E,uation (18) ensures that r^ and r^ onl, assume values within 
the Interval (0.1). For the cate where the minimum concentrations are sero and 
the maximum concentrations are th. initial values, the combination of equations 

(14) and (18) yields: 


where Cp 


» k (t) (rpCOp) VqCjJ) 

Cq = initial concentrations of fuel and oxidizer, respectively 
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2.3 Mean Turbulent Arrhenius Reaction Rate Constant and 
Mean Turbulent Reaction Rate Expressions 

In a turbulent reacting flow, the mean turbulent Arrhenius reaction rate 

constant can be calculated by treating the temperature as a continuous random 

variable and specifying an appropriate probability density function for the 

temperature. This concept is employed in Reference [4lf along with the 

analogs of equations (4) and (17) in the present work, to yield: 

kTtT = 0/ * k(t)p(t)dt (20) 

where k (t) * mean turbulent Arrhenius reaction rate constant 
The above expression provides a direct, relatively simple method of taking into 
account the effects of temperature fluctuations on the Arrhenius reaction rate 
constant in a turbulent, reacting flow. , 

Similarly, the expression for the mean turbulent reaction rate can be 
developed from an extension and combination of equations (10) and (19). Thus: 

Wp 0*^^ Wpp(t,rp,rQ) dt drp dr^ (21) 

where w- = fiaean turbulent reaction rate of fuel 
F 

This expression utilizes a joint pdf for temperature and species. Equations (20) 
and (21) can be applied to one-step mechanisms or, multi-step mechanisms, by 
considering each elementary reaction separately. The pdf*s used here are 
considered to be valid at an instant in time, and thus are not functions of time. 


ORlGtNAl. 
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In order to compare the ifiagnitude of a mean turbulent Arrhenius reaction 
rate constant to the corresponding laminar term, an "amplification ratio" is defined. 
The laminar Arrhenius reaction rate constant is calculated by inserting the mean 
dimensionless temperature in equation (17). The reaction rate constant amplifi- 
cation ratio is obtained by dividing equations (20) by (17): 


Z = 


( 22 ) 


where Z = reaction rate constant amplification ratio 
= mean turbulent reaction rate constant 
■: corresponding laminar reaction rate constant 
A similar term may be defined for the ratio of a mean turbulent reaction rate 
to the corresponding laminar reaction rate. The laminar reaction rate is 
calculated by inserting mean values of dimensionless temperature and concentrations 
into equation (19). The ratio of equation (21) and (19) defines the reaction rate 
amplification ratio: 


where Z' ® reaction rate amplification ratio 

.8- s mean turbulent reaction rate 
w 

t 

w^ B corresponding laminar reaction rate 
Equation (23) is the reaction rate amplification ratio in which the mean 
turbulent reaction rate is calculated with consideration of the combined effects of 


1 


] 



1 


12 
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temperature and species concentrations fluctuations. In order to compare this 
amplification ratio to one which considers only temperature fluctuations, the 
following is defined with the aid of equation (19): 

^ ^ yrpCV^roCy 

k(?pCp"(?oCj,)'" 

This expression involves a mean turbulent reaction rate in which a mean turbulent 
reaction rate constant is employed along with laminar values of concentrations. 

The corresponding laminar reaction rate also contains these same concentration 
terms. Thus, the concentration terms cancel. The reaction rate amplification 
ratio for this case is: 

2' = k/kj (24) 

This amplification ratio may be compared to that calculated by use of equation (23). 
Hence, combined effects of temperature and species concentrations fluctuations 
on the reaction rate may be compared to effects of only temperature fluctuations. 
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3.' Joint Probability Density Functions For 
Temperature & Species 


3.1 Most-Likely Bivariate pdf for Two Species 

A model which accounts for the combined effects of temperature fluctuations 
and concentrations fluctuations of both fuel and oxidizer species» explicitly, 
on the mean turbulent reaction rate is presented in this section. As an initial 
step in the formulation of such a model, temperature is treated as an independent 
random variable, and as such, has no effect on the species concentrations. 

Mathematically, the three-variable pdf for temperature and species is 
expressed as: 

p(t,rp,rQ) = f(t)g(rp,rQ), 0< t < 1 ( 

0 i rp 1 

0 < Tq ^ 1 

where f(t) = a pdf for temperature 

g(r_,r^) = a joint pdf for the concentrations of fuel and 
^ ® oxidizer species 


for the case where temperature is treated as an independent continuous random 
variable. Equation (25) is a valid pdf since it satisfies the following extension 
of equations (5a) and (5b) for a three-variable pdf: 

p(t,rp,rQ) i 0 

P<t.rp.rQ)dl dry dr^, = 1 
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The most-likely bivariate pdf is used as the joint pdf for the concentrations 
of fuel and oxidizer species in this model. The most likely bivariate pdf of 
Reference [41 can be written in the generalized form: 


^ n 

P (tp ) » q (tff '> exp (^1^ ) 

where q(’^ ) is the a priori probability which, for non-reactive scalars, is a 

constant. The ( X +D coefficients A are uniquely determined from the 

X 

moments and from the condition that the pdf integrates to unity. When the 
first three moments are known for the pdf, it shows excellent agreement with 
experimental data. 

As an initial step in the utilization of the most-likely bivariate pdf for 
two species, a pdf based on the first moment about the origin of each species 
concentration and the covariance of the two concentrations is selected. The 
expression for the chosen form of the most-likely bivariate pdf for two species 


(27) 


is: 


gCrptr©) = q • exp ( ^ 2 ’‘o ^ 3 ( 

where q * a constant 

r~ = dimensionless concentration of the fuel 
F 

r^ s- dimensionless concentration of oxidizer 
The continuous random variables in equation (28) are treated as passive scalars. 
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With this Simplltication, the "q" term in equation (28) is a constant. The value 

of this constant "q" term is shown in Appendix A of Ref. 2 no have no effect on the 

calculated value of the most-likely bivariate pdf and the values of the resulting 

reaction rates. 

The values of the constant coefficeints, Xq, X j,X 2* ^3* *** equation 
(28) are obtained from the simultaneous solution of the following constraint 
equations for known values of rp, r^ and r^r^: 


o/'o/' 

g(rp,rQ)drpdrQ » 1 

(29) 

o/'o/' 

rp g(rp,rQ)drpdrQ = rp 

(30) 

,j\f' 


(31) 


(rp-Fp) (ro-FQ)g(rp,rQ)drpdrQ = r^r^ 

(32) 


The pdf selected for temperature is either the beta pdf, the one variable 
most-likely pdf. as presented in reference 12] or the ramp pdf. based on certain 
criteria specified in Reference [5]. The development and utilization of these 
criteria are discussed in Reference I5l and will not be discussed here. 

Equation (25) is used in an expression for the mean turbulent reaction rate 
by considering the one-step, irreversible reaction, as in Section 2.2: 

F + V O -*> P 

where F.O.P = fuel , oxidizer and product species, respectively 

X. s stoichiometric coefficient (molar.) 
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Assuming that both n and m are equal to unity in equation (19), the reaction rate 
of the fuel may be written as: 

Wp « k(t)(rpC®p) (r^C^) (34) 

For the case of statistical independence between temperature and species 
concentration, the expression for the mean turbulent reaction rate of the fuel 
species is obtained from a combination of equations (17), (21) and (25)> This 
yields: 

"w^ a -A o/ (kjt+k 2 )®exp [“T^/(kjt+k 2 ) ] f(t) df (35) 

The corresponding value of the laminar reaction rate if determined by inserting 
the appropriate values ofT, Fp and Fq into the combination of equations (17) 
and (33). This yields: 

w^ = -A (kjT+k2)®exp [-T^/(kjT + k2) ] FpCp r^C^ (36) 

The reaction rate amplification ratio is obtained by dividi equation (3S) by (36). 
This yields: 

(kjt+k2)®exp t-T^/(kjt4k2) ] f (t) dt 
(kjT+k2)^exp t-Ty^/(kjt+k2) ] rp 


Z' 


(37) 
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Th. .mplificatlon ratio in equation (37) may be eapreased as the product of a 
term which accounts tor the effects of temperature fluctuations. ZJ . and a term 
which accounts for the effects of species concentrations fluctuations, z; . 
These terms are expressed as follows: 


exp [-T^/(k^uk 2 ) 3 f (t) 

(kjt+k2rexp [-T^/(kiT+k2) 1 




TpFQ gCrp.rp) dfp 


(38) 




Comparison of equations (3s) and (22) illustrates that z; is equal to Z. 

In this study, values of dimensionless mean concentrations are varied 
over the range 0 . 0 to 1 . 0 . and the values of mean square fluctuations of 

dimensionless concentration are varied over the range o.O <» 0.1 • The 
value of the correlation coefficient between the fuel and oxidiser is assumed to 

be -0.9 111. The covariance can then be determined from; 


■"f To ■ ’’ ’’o 

Where p * correlation coefficient 
^ (nean square fluctuation 


(39) 


The com 


putational procedure to determine is described below: 
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1. Select values for Tp, Tqi Tp and Tq . 

2. Calculate the covariance, using equation (39) and a correlation 
coefficient of -0.9. 

3. Use Newton's Method to solve for the constant coefficients of 
equation (28). 

4. Integrate equation (38) to calculate . 

The entries to the augmented matrix for Newton's Method and the related 
mathematical limitations for the most-likely bivariate pdf will be described in the next 

status report. The data obtained from this parameter study are presented in 
their entirety in Reference [^]. Representative results are presented in 
Section 4.1 . The results of this parameter study are being compared to 
experimental data, to judge the validity of this model and the assumptions on 
which it is based. 

The model given by Z = Z| * Z^ is currently being examined in a 
large-scale computer program at NASA Langley. The portion of the reaction 
rate amplification ratio which accounts for temperature fluctuations, is 
calculated from either the beta pdf or the ramp pdf . The selection criteria is 
specified in Reference 151. Z^ , the portion of the reation rate amplification 
ratio which accounts for the effects of species concentration fluctuations is 
supplied in the form of a subroutine for the computer program . 

Subroutine F1ND2 is called to find the value of the unmixedness multiplier, 

2^ given the mass fractions of species A and B and the corresponding 

fluctuations of the mass fractions. 
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The values of are stored as elements of the matrix ABCD. For 
given values of RA, RB, RSA, and RS . find Z{.. The problem is symmetric 
in values for RA and RB ; the same is true for values of RSA and RSE. 

A finite number of values of Zj. are stored in ABCD. These correspond 
to RA = 0, .1, .2, RSB = 0., .01, .02, ..,.1. 

The given values of RA, RB , RSA, and RSB are decimal figures with 
more than one significant figure. Therefore, interpolation is necessary. This 
is accomplished by finding the l6 values of Z^ which are nearest neighbors to 
the desired value - one larger and one smaller. For each given value of 
RA, RB, RSA, and RSB, these sixteen values are added and the result is 

divided by sixteen. 


3.2 An Alternative Three— Variable Model 

A full three variable model for temperature and two species concentrations 
is discussed in this section. The assumption of statistical independence between 
temperature and species concentration has been dropped in this model. The 
assumed pdf for the full three-variable model uses a joint pdf for bo^ the 
temperature and species concentrations, and is of the form: 



where q = a priori probability 

X^-Xg = adjustable constants 
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This full three-variable model is an extension of the statistically •’most- 
likely” bivariate pdf in Reference l4l as shown in equation (27). when all of 
the first and second moments are known for the pdf . The physical justification 


for using this pdf and extensions of it was previously discussed in Section 3.1 . 

The values of the constant coefficients, ^2* 5* ^6’ 


in equation 
constraint 


(40) are obtained from the simultaneous solution of the following 
equations for known values of t, *“p» ^o* 


o/‘t/‘o;‘rp(t, rp,rQ)dtdrpdrQ - Fp 

111 _ 

0/ 0/ of r^pCt.rp.rQ) dtdrpdrQ * Tq 

o/S/^o/* (t -T) (rp - Fp) p(t,rp,rQ)dtdrpdrQ = fr^ 

0 / \ / 0 ; ' (t -"t) (rQ - Fq) p(t . rp . r^) dtdrpdr^ = T^ 

o/'o/'o/'(rp-Fp)(ro-FQ)p(t.rp.rQ)dtdrpdrQ = 7 ^ 


(41) 

(42) 

(43) 

(44) 

(45) 

(46) 

(47) 


The temperature and species concentration in equation (40) are treated 
as passive scalars. With this simplification, the “a priori” probability q is a 
constant . The proof shown in Ref. 2 can be extneded to show that the q 
term in equation (40) will not effect the calculated value of the pdf, p(t,rp,rQ), 
and therefore has no effect on the resulting reaction rate. 
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The pdf, p(t,rp,rQ), is used in an expression for the mean turbulent reaction 
rate by considering any irreversible reaction of the form: 

A + B C + D ' ^ 

where A, B = reactants species 
C, D = product species 

Assuming that both n and m are equal to unity in equation (14), the reaction rate 
of the fuel may be written as: 

ip . k(t) (rpCp (roCj») 

A combination of equations (17). (21), and (25) yield an expression tor the 
mean turbulent reaction rate as: 

w^ = - A CpCou/%/'o/ ^ 

• rprQp(t,rp,rQ)dtdrpdrQ 

The corresponding value of the laminar reaction rate is determined by inserting 
the appropriate values of Fp, Fq. andT into the combination of equations (17). 
and (33). This yields: 


w, * -A(k,t+kJ exp ( 


-T . /(k,t+k«) ) 


2 ' ' F F O O 


The reaction rate amplification ratio is obtained by dividing equation (50) by 


(51). yielding: 


i. 
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o/*p/*o/'rprQ (kjt+k2)®exp [-T^/(klt+k2)] P (t, rp.rQMtdrpdr^ 


(52) 


F^Fq (k^T+k2)®exp[-T^/(k jt+k2)] 


In this parameter study, the dimensionless mean temperature and the 
dimensionless mean concentrations are varied over the range from 0.1 to 0.9. 

The values of mean square fluctuations of dimensionless concentrations and 
mean square fluctuations of dimensionless temperature are varied over the range 
from 0.01 to 0.09. The corresponding value of the covariance for r^r^ 
is determined from equation (39) using an assumed value of the correlation 
coefficient of -0.9^^^ The covariance for and r^f are also calculated 
with a correlation coefficient of -0.9 for a similar form of equation (39). 

The computational procedure to determine the values of the reaction 
rate amplification ratio is described below: 

Procedure 

- - - — 7 ? 2 " 

1. Select values for t, t*p, ^ 

2. Calculate the three covariances using equation (39) and a value of 
the correlation coefficient of -0.9. 

3. Use Newton's Method to solve for the constant coefficients of 
equation (40). 

4. Choose values for B, T^^, 

5. Numerically integrate equation (52) to calculate Z . 

The values chosen for A, B, T^, C* , C* are given in Table 1. 

Selected results of the parameter study undertaken are 
presented in the next section. 
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A> Results 

4.1 Selected Results of the Most-Likely Bivariate Pdf For Tvo Species 

All numerical integrations are performed with the use of Simpson's . 

Rule. The number of subintervals on the integration interval of [0,1] are large 
enough so that no truncation error is incurred. The resulting value of Z^> 
is accurate to five places. 

In addition to the limitations of the program from Newton's Method, 
convergence problems were found from three other causes: 

1 . The value of ( X2 + ^3 was a term close to zero in the denominator. 

OQ 

2. The value of the pdf became greater than 1 x 10'^ , or less than 

•^8 

1 X lO'^ . This can be caused from a poor initial guess, or 
extreme values of the parameters. 

3. The last problem was from one of the constraint equations. Equation 
(32) can re rewritten as: 

« / ' 0 / ‘ p'''f ’ ■'o’ ' ■'f ''6 * Vo 

Since rj-r^ is always a negative term, for certain values of the parameters 
the right hand side is less than zero. This violates statistical theory and 
hence no convergence is possible for these values. 

Figures 1-4 show the variation of 2^' with increasing mean square 
fluctuations at constant mean species concentrations. As the fluctuations increase, 
Zr’ always decreases. The lower limit ofZ r' is zero, the case where there is 
no reaction. This will occur when the turbulence is so great that the fluctuations 
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inhibit the reaction, or where one of the species is not present. 

Because of the convergence problems indicated with the computer program, 
results for high values of the mean square fluctuations were linearly extrapolated. 
Once a value of zero was reached , it remained zero for all higher values of mean 
square fluctuations. For the upper values of mean species concentrations Zr' 
is approximately constant and approaches a value of 1.0. 

Figure 5 shows the variation of Zr' vs mean species concentration for 
constant values of mean square fluctuations. The value of Zr' increases with 
increasing mean species concentration. 


4.2 Preliminary Results for the Full Three-Variable Model 

To decrease the computer time to calculate each Z value for the full 
three-variable model, an exact transformation was used. The triple integrals, 
which must be calculated to use Newton’s Method for non-linear systems, where 
reduced to double integrals where the variables are separable. All of the double 


integrals have the following form: 

^ -a/4)o 

f f(T)/ o^e 

^ o=f(T) 

where N s 0 , 1 , 2, 3> 4 


dT 


( 54 ) 


The augmented matrix for Newton’s Method and a sample derivation of a double 
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integral will be in n future progress report. When N O, the integral over 
o is a scaled Version of the bell shaped curve. For N «a 1, 2, 3, or 4. 
the integral over o is a function of o and tne bell shaped curve. 

One of the mathematical difficulties of Newton's Method is that convergence 
of this method is dependent on an initial guess within acceptable bounds of .the 
exact answer. 

As a first step toward an initial guess subroutine, restrictions wert 
placed on the variables to reduce the full three-variable model to the two-variable 
model. (Note: In the following discussion, t is dimensionless temperature 
and r and x are dimensionless mean species concentrations A and b). 

-.2 

1 he restrictions on the three-variable model arcT - r and f r’ . 

t and r can bo considered dummy variables since they are integrated over tl,e 
same limits. Therefore let t r y. The pdf can now be written as: 

2 

pdf cxp(j ^ 

• expd ^ ^ I ^ 5 ^ 

This rrow gives a value for the pdf which can be evaluated from the one and two- 
variable pdfs. The initial guess subroutine can be developed by using powei 
scries Around exact solutions. 

I„ ll„. wlM-r.- Z - Z, Z^. t is not o luiu lio.i of either r ef x ls|>.'.'ies 

A Of 10. Theretoff, the inteRial over t con be pineeil inside the double Intejtral 


t.'e’.e 
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over r and x and Z can be calculated from: 


Z^of\f\f'Q exp(Ag + xyt+X^’ 

• expCxjy + t + X 2 X + X jTx) dxdrdt (56) 

where Q * rx exp [ 3 + ^ 2 ^^ 

1 2 

Z for the full three variable model under the restrictions indicated for 
equation ( 55 ) can be calculated from: 

I 1 2 

z = 0/ 0/ exp [| Xq +x 1 y + X4 y ] 

•exp fix 0 + X2 y + X3 X + ( X 5 + X 5 ) yx 1 dydx (57) 


Comparison of the Z value for the three variable model and Z calculated 
by multiplying the one-variable temperature only pdf by the two-variable bivariate 
pdf for two species can be evaluated by taking the difference between the values 
and determine if it can be made arbitrarily small. In equation form: 

, ,, . . L Xo ♦ Xjl 4 ^rx 

Z - =o/o/o/\e 

- expC X 3 + X i't +X 2 ^ Xj)+ Xl r 4x ^ x +x grx:^dxdrdt 

(5S) 


Under the restrictions of t = r 


and 



the difference between 


the Z values is exactly zero, if the correlation coefficients are as follows: 



P 


rx 


+ 1.0 

-0.9 

-0.9 
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5. Conclusions and Direction for Future Work 

5.1 Conclusions 

5.1.1 Effects of Species Fluctuations on the Reaction Rate 

The effects of species concentration fluctuations on the reaction rate are 
assessed by treating the species concentrations as continuous random variables. 
A joint pdf is used to relate the variables. The total amplification ratio is the 
mean turbulent reaction rate divided by the corresponding laminar value. The 
reaction rate amplification ratio, which only accounts for the effects of species 
concentrations, is discussed here. The results obtained for a parameter study 
on zf are the following: 

1 . The value of Zr is between 0 and 1 .0. 

2. For constant values of mean species concentrations of A and B 

» ?«) and for constant fluctuations of species A (?' ), Zr 
A o — ^ ^ 

always decreases as fluctuations of species B(Fg ) increases. 

3. For constant values of fluctuations for species A and B and _for 
constant mean species concentration of A, Zr increases as r^ 
increases. 

5.2 Direction for Future Work 

1. Having exact answers to the three variable model, evaluate the 
error in the program and write an initial guess subroutine. 

2. See how Z compares to Zt * Zr for all of the mean quantities 
not equal. 


3. Compare Z vs Z^- for different correlation coefficients. One 

suggestion is te set the value of the correlation cwfficlent between 
species equal to -0.9 and between temperature and species equal to 
-0.99. These values were suggested by Antaki [11. 
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ihe ttudy of hydrogen fuelea super&onic combustion 
(scr^-'^'et?) has been a major project in the research programs 
of the Isnfley Research Center. A large computer profra-^ was 
developed which solves flowing or static chemical kinetic prob- 
le-'S involving many chemical species (Kef. 7). Ihe program uses 
the conservation equations (Ref. 8) for one-aimensional steauy- 
state flow -.o analyze a complex reacting gas mixture flowing 
through an arbitrarily assignea area. 

It is known that the computational time requirem.ent for 
this com.ruter program, is proportional to the num.ber of species 
^nd reactions beinr treatea. it is the purpose of this task • 
to develor a methoa which reduces these num;uers in the stua;,- of 
hydrerer..'- ir co'-ousticr. and at the same timie preserve the correct 
r'-vsic'^ l-che’’ioa 1 oeh-.vior. .he result of this »rfork will be to 

reduce computer time requirements. 

.he method that Vtill oe usee in this p:.per will ce -n 
c ' on and refininr of the metr.oa protoseu oy w.ninitr / 

•'ethod involves the tracking of a "trigter" specicc v.nic:. 
ic usee to determine whether the flo'w is in an " ignition" m.o..v 
or in - "combustion" moue. ignition ana com.Dustion are aefi.-.v.. 
in the classical sense. Ignition delay tim.e is t' ken to it t:.t 
time required for the temperature increase to reach five rercer.a 
of the overall temperature increase. 


‘ if-nitinr. 


nn 


itial ‘equilibrium ~ ‘initi: !^ 
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.Vhen the teT.perature is less than flov/ is in the 

••irnition" tpode. otherwise, the flow is in the •• combustion" 
mode (see firure 6). 

This work '*^0 while an extensive 

chen'ical packare is needed to aescrioe the "ignition" mode, a 
s-aller TDackage is sufficient to deal with the "combustion" mooe. 
iV'r.ns -nd Schexnr.yaer concluded that a reaction system involv- 
inr 12 species (referred to as ^5(1^^)) w-s required to oescrioe 
"irnition" processes while an fa(7) system sufficed to deal with 
the "com.bustion" mode (Kef.lO). This paper will try to minimize 
the co"'t'Ut?-tiona 1 time requirements once the "comcustion" process 
is initiated (i.e. a sm.aller package or one that takes less compu- 
tational ti-"e than the c(7) system might suffice to deal witr. 
+VC "comcustion" moue). 

1 ^7(1 3) system will be used as the "test" chemical kinetic 
r-c>a-"e (lable 2). Any system with sm.aller numoers will be 
"tested" ar^inst the 37(13) system to see if the correct rhysic^-l- 
cwe"ic?l behavior is preserved, dhis will oe a cco-^rliEneu in tr.e 
followinr way First, the full 37(o) system will oe run fro- the 
initial state to equiliorium. Next the test system, will oe run 
fro- the point where the "com.bustion" in the 37(13) system initi- 
ates. xt.e test system.'s initial values (pressure, tem.t er- ture , 
equivalence ratio, chemical com.position) will oe that of tr.e 
37(13) system's values at the point of combustion. Temperature- 
time profiles ano chemical behaviors will be compared (see figure 
7 ). If the difference between the 37(13) system, and the test 
svste- is within an acceptaole range then the test systor^ will uc 
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sufficient to der. 1 with the ••combustion" process. 

An 8(7) system (Vable 3) will be tested to confirm the 
results of Evans and Schexnayder, In oraer to further reduce 
comput- t i - 1 time, a ^(d) global model ('.'Lie 4) developed by 
Chinif^ and Ro, *rr (hef.ll) will be tested. Lastly the p^per 
will try to test the partial equilibrium assumption (Ref. 12 ) 
on the t oove systems in an effort to absolutely mirii"ire compu- 
tational time. .lO accom.plish this, a large num.ber of one-dime.n- 
sion'l constant pressure, - air computations will be perfor.m.e 
in the ranges 

C.5 ^ 1.5 

{•5C 4 r 4 IcOQ 

0.5 ^ F - i-C 

where / is the eq'uivulence ratio, 'i the temperature in , ^nu 
r t' e pressure in atmospheres, a list of tnese inciviaual cases 
is riven in 'rable 5 . 

1 lots of the mole fraction of the trirrer species vers'ur 
irnition time and ignition time versus ignition tempers ture v.ill 
re ’^ade in accordance with the process described in reference ?. 
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Frelir,in&ry hesults 


The ’Tiajor prohlen. in the first serr:ester of work he.d 
been in the orerf.tion of the computer progran'. The progr^irr 
had been runninr smoothly on the N. Y. 1. operating system 
until the system was switcheu to a <Jior,R. At present the 
r»-o r"'f is not working out will oe operating within the 
nf>" As 3 result progress has oeen slO'- . 

.'■<} first 35 of the 45 cases have been run using the 


extensive 37(13) package. Cases io-21 using this package 
were terminated cue to excessive computer running times. 


rrcblem h- s been encounterea before in other stucies. 
?i:^teen cases have ueer. run using the b(7) fro- ir.ni 
‘ libriu-. Figures fi- 12 s.how the o(7) syste'r com.p: reu 


icr. 

to 


the ”'■( 13 ) system, eta rressure of one atmosphere, an equiva- 
lence r' tie cf 0.3, ar.u various te.mperatures. ^nis set cf 


co-r' 


: SO! 


is tvviC' l 0^ the results founu in co.mp;- ring the 


s'‘stem v.'ith tr.e g7(^3) system for the entire prouesc.. .■% 
■ :nc '^OO^K the. 8(7) syste.m; lags oer.ina the 37(13) systt'- 


-A'^rera tures of i0''0°h, llOO^h, anci lii00°n the ?7(13) 

lags behincj t!.e c(7) system., ih.erefore, oetv.een r. 

ana a transition tikes place. As a result there is a 

te'-rera ture range at approximately 950°K where the- t(7) sycter 
descriser doth the ignition ana com.bustion processes. 

?i--ure 13 shows the comparison of the 37(13) ana the 
•^”) svr‘e-s in the coruusticn moce. As shown by bvans anu 
uC' f xir yaer the b(7) syste- uescrioes the com.iustion procer:.- 
'ore runs of the b(7) system starting at the oom.uustion pcii.t 
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r^f the 3?(13) sycter are neeuea for further verification. 

Fic'ureii}, shows a logarithmic plot of mole fraction o 
G'i versus ignition time. As shown in reference 9 the plot 
is linear and is thereby a g.ooa trigger species. Kore tea 
cases will be needea to verify this. 

Figure I 5 shows a semilog.arithrr.ic plot cf ignition 
tir-e versus the reciprocal of initial temperature for 
both, the &(7} and ?7(i3) systems. Ihe results here confir 
the work of Rogers and Schexnayaer. 
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Conclusions ana coservtt tions 


Vhese are preliminary conclusions basea only on the 
initi^?! vs'ork done so farj 


1 . 




•'he ^’U7) system does not describe the entire process 
(i.e. irnition and combustion) except at 950®h, one 
? tmoschere , and an equivalence ratio of 0 . 5 - 
"he b ( "^ ) system oescrioes the comoustion process, 
it one atmosphere ^nd an equivalence ratio Oj. 0.5 
the lo£^'rith.'"ic plot of mole fraction of Or. versus 
irnition time is linear. 


i.t the sa^e conditions as aoove, the 
temnera ture on irnition time for the 
confirms the result- of reference S- 


effect of initic 1 

3 ?( 13 ) pacKace 


37 


FLAME STABILIZATION STUDIES USING A 


PERFECTLY-STIRRED REACTOR MODEL 


ORIGirvAl. PAQE IS 

OF POOR QUAUTY 


The use of the stirred reactor concept to model combustion 
regions of intense mixing and recirculation dates to Longwell, 
et. al. (Refs. 12 and 13). Longwell's laboratory stirred 
reactor, along with regions in a SCRAMJET engine which may be 
modeled using stirred reactor theory, is shown in Pig. 16. 

Reactors in which reactants and inerts exist in the gas phase 
only are termed "perfectly stirred"; multi-phase systems are 

"well-stirred." 

A version of the perfectly-stirred reactor code developed 
by Pratt (Ref. 14) was used along with chemical kinetic mechanisms 
developed at LRC for hydrogen, hydrocarbon and silane combustion 
(e.g. Refs 8, 15, and 16) to investigate blowout limits for neat 
fuels and for mixtures of these fuels. Results are shown in 
Figs. 17 through 27. The blowout limit, in all cases, was deter- 
mined by examining computed curves of reactor temperature as a 
function of the mas® flow rate per unit reactor volume, in/V. The 
blowout limit was then determined as suggested in the sketch below 

T ^ 



Equilibrium Flame 
Temperature 


Specified Reactants, 
Initial Temperature, 

System Pressure, and 
Fuel-air Equivalence Ratio 
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The residence time In the renctor, related to ir»/V by 

tp * p/(tu/V) 

is simultaneously calculated in the computer code. (Note that a 
recent criterion for blowout developed by R. C. Rogers was not 
available at the time the calculations reported upon herein were 
made . ) 

Fig. 17 illustrates the effect of initial reactant temperature 
and quantity of silane in the fuel on minimum residence times 
required to avoid blowout of the flame. These results are for 
H 2 and as the principal fuel constituents. It is of interest 
that increasing SiH^ in the fuel results in a more substantial 
effect for H 2 than for CH^. On the other hand, the results for H 2 
suggest that, at high initial temperatures, little stabilizing 
effect is achieved for H.,/SiH^ until the silane concentration exceeds 
5% by volume. Increasing the initial temperature, as expected, 
increases the stable flame region. These same results may be noted 
by reference to the cross-plot in Fig. 18. 

Figs. 19 and 20 illustrate the flame stabilizing effect of 
increasing the system pressure. As can be seen, the effect is sub- 
stantial going from 0.2 atm to 2.0 atm wherein, at 1000®K for 
example, minimum required residence time decreases from about 7 ;isec 
(at 0.2 atm) to about 1 iisec (at 2.0 atm). 

The effect of equivalence ratio is indicated in Figs. 21 and 22 
for a 20% SiH^/80% H 2 fuel (by volume) at 1.0 atm. The stabilizing 
effect of increasing the equivalence ratio C is seen to be predomi- 
nant at thf' lower values of initial reactant temperature T^. For 

-4 

example, at T » 600*K, t« decreases from about 3 x 10 sec (4 - 0.2 




to about 3 X lO"® sec ((|) » 2.0); a decrease of two orders-of- 

magnitude. On the other hand, at • 1200 the decrease in 

t« is from about 2 ysec (4» = 0.2) to 0.5 ysec (♦ - 2.0). 

R 

The traditional stirred reactor plot of 4> versus in/V is 
shown for SiH^/Hj fuels (T^ - 600«K, p * 1 atm) in Fig. 23. The 
stabilizing effect of increased silane concentrations is again 
noted. Of particular interest are the following; 

1. The locus of maximum £n/V values increases from (|i = 1 for 
pure H 2 to about = 2.5 for pure SiH^; and 

2 , ^ reversal occurs below ij) ® 1.0 wherein the stable flame 
region is greater in extent for lesser quantities of silane in the 
fuel . 

A possible explanation for these two findings is suggested in 
Table 6 where it is noted that the principal initiating reactions 
for silane oxidation are SiH^ pyrolysis followed by the very rapid 
oxidation of the pyrolysis produce SiH 2 » As seen, these two steps 
are exothermic. However, the direct step ©2 oxidation of SiH^ to 
SiHj and HO 2 is mildly endothermic when compared with the highly 
endothermic steps which initiate the CH^ and H 2 oxidation chains. 
Hence, high silane concentrations relative to oxygen (i.e., high 4 ; 
values) tend to result in greater flame stability due to the 
exothermicity and relative rapidity of the principal chain-initiating 

steps. 

On the other hand, below 4 * * 1 flame stability is favored by 
the higher exothermicity of the overall reaction (Table 6 , 4'<D • 

Since the heat of combustion of silane is substantially below that 
of hydrogen and even the hydrocarbons, at low equivalence ratios 
increased silane concentrations appear to have a destabilizing effect 


A similar result is seen in Fig. 24 where SIH^/CH^ blowout limits 
are shown. Fig. 25 again suggests the picture detailed above. For 
example, the activation energy of the reaction 
CH^ + M -► CH 3 + H + M 

is given in Ref.H as about 85,755 Kcal/gmole, while for 
CjHg + M CH3 + CjHg + M 

it is given as 64,968 Kcal/gmole. Hence, the relative rapidity of 
the latter tends to favor C^Hg at <|)>1. On the other hand, the 
higher heat of combvistion of methane favors it at 4x1. These 
predicted trends are verified in Fig. 25. 

Figs. 26 and 27 summarize the results obtained for the neat 
fuels in terms of residence time (Fig. 26) and ifn/V (Fig. 27) . In 
general, for 4>>1 flame stabilization is favored by the use of 
silane in the fuel mixture and by the use of the higher hydrocarbons 
(above CK^) . For 4><1, silane may serve to reduce flame stability 
and the use of methane is favored over the higher hydrocarbons. Also 
favoring flame stability are high initial reactant temperatures and 
high pressures. 
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. j ♦ 


( 


C 


A 

B 

C 

C® 

E(x) 

P 


y) 

y) 


gTxTyT 


h(x) 

k 

ki , k2 



m 


max 

min 

n 

0 

P 

p(x) 
pCx, y) 
r 


pre-exponential constant 
temperature exponent 
instantaneous concentration 
initial concentration 
expected (mean) value of x 
denotes a fuel species 

joint probability density function of (x, y). 

function of the two-dimensional continuous random 
variable (x, y) 

mean value of g(x, y) 

function of the continuous random variable, x 

Arrhenius reaction rate constant 

constants dependent upon T . and T 

min max 

laminar Arrhenius reaction rate constant 

mean turbulent Arrhenius reaction rate constant 

constant, defined in equation ( 30 ) 

denotes maximum value 

denotes minimum value 

order of reaction constant 

denotes oxidizer species 

denotes a product species 

probability density function of x 

joint probability density function of (x, y) 

dimensionless concentration and constant defined 
by equation ( 30 ) 
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r mean dimensionless concentration 

mean square fluctuation of dimensionless 
concentration 

t dimensionless temperature 

t mean dimensionless temperature 

t^ mean square fluctuation of dimensionless 

temperature 

T temperature 

activation temperature 

V(x) variance of x 

reaction rate 

laminar reaction rate 

mean turbulent reaction rate 

Z reaction rate constant amplification ratio 

Z' reaction rate amplification ratio 

r gamma function 

u mean of a continuous random variable, x 

X 

u, kth moment of a random variable 

k 

o2 variance of a continuous random variable 

p , correlation coefficient 
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Figure 4: Variation of the reaction rate amplification 

ratio with mean square fluctuations for 



Figure 5: Variation of the reaction rate amplification 
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